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CHAPTER 1 


INTRODUCTION 

It is clear from previous work [1-4] over the past six years that most, if not all, 
of the lightning strikes encountered by the NASA F106B thunderstorm research aircraft 
were triggered by the aircraft. This places increased importance on the investigation 
of the triggered lightning environment. The central question to be answered may be 
stated very simply: "Under what conditions will an aircraft trigger a lightning strike?" 
The answer to this question depends on a number of factors, the most important of 
which are listed below. 

1) Electric field and charge distribution within the thundercloud, 
particularly near the aircraft; 

2) Shape of the aircraft; 

3) Orientation of the aircraft with respect to the ambient electric field; 

4) Altitude of the aircraft; 

5) Net electric charge on the aircraft; 

6) Thunderstorm particle environment near the aircraft. 

The local electric field in the vicinity of the aircraft is most influential in the 
generation of a triggered strike. It must be above a certain threshold in order for 
electrical corona to appear initially, and after a leader has developed and is 
propagating away from the aircraft, it serves to maintain that leader. The charge 
distribution within the cloud is, of course, responsible for the electric field distribution, 
and also provides a source of charge for K changes, or intracloud return strokes. 

The shape of the aircraft is instrumental in enhancing an ambient electric 
field which is below air breakdown intensity up to a level at which air breakdown can 
occur. This normally occurs at locations having large field enhancements, such as the 
nose boom, wing tips, and tail of the FI 06B. 
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Orientation of the aircraft with respect to the ambient field is important 
because field enhancements are, in general, different for different directions of the 
field. As an example consider the F106B in level flight with an ambient electric field 
oriented from wing to wing. The largest local fields will be located at the wing tips of 
the aircraft, and may be below air breakdown intensity. If the aircraft now turns so that 
the ambient field is oriented along the fuselage, the largest local fields will be located 
at the nose. Because the field enhancement for this latter orientation is larger than for 
the former, it is possible that a triggered strike could be initiated. Hence some 
orientations of aircraft and field may be favored in triggering a lightning strike. 

The altitude at which an aircraft flies influences the triggering process in that 
the minimum air breakdown field is a function of air density. The functional 
relationship is roughly linear with the breakdown field decreasing as the air density 
decreases. Therefore, for a given ambient electric field, triggering of lightning 
becomes more likely as an aircraft climbs in altitude. 

In normal flight conditions an aircraft acquires a net electric charge through 
collisions with particles in the air. This net charge, if large enough, can significantly 
alter the ambient field level at which triggered lightning can occur. The alteration can 
be in either direction. Depending on a number of factors, one sign of charge may 
lower the triggering field, and the opposite sign of charge increase it. 

Thunderstorm particles, in addition to providing a source of net charge on an 
aircraft through collisions, generate local volumes of enhanced fields. These particles, 
particularly ice particles, are irregularly shaped and conductive on static time scales. 
They can produce field enhancements of ten or more in a small volume surrounding 
the particle. This field enhancement in turn produces local corona. A large number 
density of these particles may then significantly alter the effective air breakdown 
characteristics in a region of space. 

This report investigates various aspects of the triggered lightning 
environment. In Chapter 2 an effort is made to measure the form factors of the F106B 
aircraft with the use of a metallized plastic model placed in a uniform electric field. The 
measured form factors are compared with analytic and numerical form factors to 
validate the measurement technique. The form factors are important in that they allow 
one to determine the minimum triggering electric field for an aircraft. They may also be 
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used in the analysis of electric field mill records to derive ambient electric fields from 
local fields measured on the aircraft. 

In Chapter 3 a two dimensional nonlinear lightning channel model is 
developed to investigate the late time (several microseconds) behavior of a triggered 
lightning event. Results from a simulated event are presented. 

Chapter 4 details the analysis of internal wires on the F106B. Responses of 
both fuselage and wing wires have been recorded during triggered events. These 
wires are modeled using a linear triggered lightning model to obtain sources for the 
wires. The wires themselves are modeled as transmission lines with discrete sources 
and 50 ohm terminations. Results compare favorably to measured data. 

Chapter 5 contains analysis relating to the effect thunderstorm particles have 
on a triggered lightning event. A general formalism is developed that allows one to 
vary size, shape, and density of these particles. A simulation is presented that shows 
the effect particles have on the responses measured on the F106B. 

Chapter 6 presents a summary of the work documented here and a recap of 
the major conclusions. 
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CHAPTER 2 


MEASUREMENT OF THE ELECTRIC FIELD ENHANCEMENT FACTORS 

ON AN AIRCRAFT 


2.1 Introduction 

The work discussed here stems from an attempt to determine experimentally 
the electric field form factors associated with locations of the field mills on the NASA 
F106B research aircraft. The experimental method chosen utilizes a small scale 
model aircraft, coated with a conductive paint, and immersed in the field produced 
between two large parallel plates. A small conductive probe is touched to the model at 
the appropriate location and the charge carried away from the surface is carefully 
measured. This is compared to the charge transferred to the probe after touching one 
of the parallel plates to get a value for the field enhancement. The method works 
because the net charge transferred to the probe from the test object is linearly 
dependent on the local electric field strength. This method has been utilized by others 
[5] and [6], but as far as is known, the accuracy of the method has never been reported. 

Prior to making measurements on the model aircraft, tests were made with 
several simple shapes such as spheres, cylinders and cubes. The sphere and the 
infinite cylinder lying normal to the field direction can be solved analytically and thus 
serve as exact checks on the test technique. The cube and the finite cylinder oriented 
in the direction of the field were analyzed utilizing a very finely gridded finite difference 
method. Comparison of the form factors obtained from the analytical and finite 
difference solutions with measurements indicate that the experimental method has 
serious flaws. The errors involved increase with the magnitude of the form factors and 
can reach values of 200 to 300% for form factors of the order of 10 (easily reached 
near the extremes of a finite cylinder oriented along the field lines). As a result of this 
investigation, the experimental method is not recommended for calibration of field 
meters located on aircraft and should not be relied upon in any application where the 
local spatial derivatives of the electric field on the model are large over the dimensions 
of the probe to be used as a sensor. 

In order to measure (or more correctly to infer) the ambient electric field in 
the free atmosphere, sensing instruments must be carried aloft by some vehicle such 
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as a balloon or aircraft. However, the presence of the vehicle and/or the instrument 
itself perturbs the ambient field, thus making the measurements incorrect. While this 
effect is indeed vexing, it does not necessarily preclude estimating the ambient field 
from measurements, since the measured field must certainly be proportional to the 
ambient field. The proportionality factor (or form factor) is a function of the geometry 
and material of the vehicle and is specific for a given point on theVehicle and a given 
ambient field direction. As a result, the task of estimating the ambient field from 
measurements taken near the conducting skin of a vehicle such as an aircraft rests 
very heavily on our ability to obtain accurate estimates of the form factor at the various 
sensor locations. 

One method of obtaining the sensor location form factors is by mathematical 
modeling. The method of finite difference solutions to Maxwell’s equations has been 
utilized for the NASA F106B research aircraft to compute form factors at the locations 
of the field mills. However, there is still some question as to the accuracy of these 
results since the computed values depend to some extent upon the cell size utilized in 
the model. Electric fields computed by this technique are averages over the 
dimensions of the cell. If fields are changing rapidly across the cell then local values 
may be significantly different than the computed value for the cell. For this reason, 
there has been a desire to check the results acquired through modeling by some other 
method. 


Another method that has been utilized in the past to estimate aircraft form 
factors utilizes a scale-model aircraft in a parallel plate electric field. The model, which 
is coated with conductive paint, is touched at the points of concern with a small 
conducting probe. The charge transferred to the probe is proportional to the electric 
field in the region where the probe makes contact with the aircraft. By touching the 
aircraft and measuring the charge, and subsequently touching one of the end plates 
and measuring the charge, one can presumably estimate the form factor by computing 
the ratio of the charges transferred. It is this method that is the subject of the 
investigation discussed here. 

2.2 The Test Configuration 

A 1 :48 scale plastic model coated with a highly conductive silver paint is 
utilized for the F106B. The eiectric field is generated between 4' x 4' aluminum plates 
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separated by exactly one meter. A potential of 5 kV is applied between the plates by a 
dry cell battery pack, thereby maintaining a balanced system to ground. An 
unbalanced system with one plate grounded through a power supply was found to 
create problems with the charge transfer to the probe and led to significantly 
unbalanced charge transfer between one plate and the other. Isolating the system 
and keeping any net charge from residing on the plates or the test object cured this 
problem. 


The test configuration initially set up for measurements on the model aircraft 
is shown in Figure 2.1 (a). Figure 2.1 (b) shows one of the test spheres mounted in the 
same parallel plate system. Figure 2.2 shows the charge measurement instrument- 
ation. A Keithley Model 616 electrometer was used for charge measurement in a 
modified "ice pail" arrangement. 

During the early measurements, it was found that a non-uniform field existed 
between the plates, giving rise to an error of a few percent. Therefore, a system of 
guard rings coupled to each other with 10 megohm resistors was installed between 
the plates and is shown in Figure 2.3. The field with guard rings in place was found to 
be uniform within a small fraction of 1 %. 

Numerical modeling of the parallel plate test arrangement was performed to 
verify the quality of the field distribution between the plates. This modeling was done 
using two dimensional cylindrical symmetry to increase spatial resolution as much as 
possible. Therefore, the plates, rather than being square as in the experimental 
situation, were round. Similarly the modeled guard rings had a circular shape rather 
than square. The area of the modeled plates was chosen to match that of the real 
square plates. A fixed potential was imposed between the modeled plates and the 
static field distribution calculated using the finite difference technique. Results of the 
analysis are shown in Figures 2.4 - 2.6. Figure 2.4 shows contour maps of the axial 
and radial electric field magnitudes for the plates without guard rings. Only one 
quadrant is shown, as symmetry exists around the axis and the center line between 
the plates. Also, the field magnitude shown is normalized to that at the exact center of 
either plate. Note in particular the distinct falloff in field near the center between the 
plates, due entirely to edge effects. 
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Figure 2.2 The Charge Measurement Instrumentation 


Figure 2.3 The Parallel Plate Configuration with Guard Rings Installed 
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Figure 2.5(b) Radial Field Mapping for Plates with Guard Rings 


12 











The effect of including guard rings with the plates is shown in Figure 2.5. 

The field distribution is not significantly changed, and may be slightly more 
nonuniform. It is necessary to connect the guard rings with large resistances to 
smooth the field distribution. This ensures that a uniform potential gradient is 
maintained at the edge of the plates. The large size of the resistors holds the actual 
current flow from the power source to a small value. The field contours for this 
situation are shown in Figure 2.6. Although the fields appear chaotic, closer 
examination reveals that the field magnitude is axial and everywhere nearly the same. 
It was therefore decided that this guard ring/plate configuration adequately provided 
the necessary uniform electric field. 

Figures 2.7(a) and 2.7(b) show the array of spheres and cylinders as well as 
the cube utilized in the method evaluation tests. The spheres range in diameter from 
1.2 to 3.3 inches while the cylinders range from .88 inches to 3.1 inches in diameter. 
The cube is 2.5 inches on a side. 

2.3. The Charge Measurement Instrumentation 

Initially, Faraday's "ice pail" method was tried as the charge measuring 
approach. The circuit is shown in Figure 2.8(a). The charge is determined by 
measuring the capacitive coupling between the charged probe and the surfaces of the 
"ice pail." However, it was found that this configuration was much too susceptible to 
external influences such as operator body movements and small charges that 
occasionally reside on the dielectric parts of the probe. Therefore, a modification was 
developed as shown in the circuit of 2.8(b). Here, the gallon paint-can (shown in 
Figure 2.1 (b)) provides a grounded shield to the sensing sphere inside the can. In this 
latter circuit, the measuring probe is inserted through the hole in the outer can and 
contact is made with the inner copper sphere. Since the input capacitance of the 
meter and inner sphere is many orders of magnitude larger than the tiny metal sphere 
of the probe, essentially all of the charge residing on the probe is transferred to the 
instrument. The Keithley 616 is operated in the "Coulomb" mode (and thus measures 
the charge transferred directly). The measurement resolution of the system as shown 
is approximately ±1 x lO' 1 ^ Coulombs, and is very stable and noise free. 

The probe configuration utilized for the measurements is shown in Figure 
2.9. The probe sensor is a 1/8" diameter ball-bearing bonded to a 1/8" diameter lucite 

15 



(a) Calibration Spheres and Cube 




(b) Test Configuration Used for the Tests 

Figure 2.8 The Charge Measuring Circuits 
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0.5 " Diameter x 3.0' Long Oak Rod 




0.125" Diameter x 6" Long 
Lucite Rod 


Tapered to 3/32" Diameter 


0.125" Diameter Metal Sphere 


Figure 2.9 The Charge Transfer Probe Configuration 
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rod. The lucite rod is tapered to 3/32" diameter at the tip where it is bonded to the 
small sphere. The handle for the probe is a 1/2" diameter oak rod approximately 3 feet 
in length. Other materials such as nylon and teflon were tested as handle material but 
proved to be no better than the wood. Apparently the lucite rod at the sensor provided 
all the electrical isolation required for the sensor. 

2.4 Measurements on Canonical Shapes 

Early in this investigation, it became clear that something was wrong either 
with the basic assumptions of the method or with the experimental technique being 
used. Initial measurement of the maximum form factor of a 1 inch diameter sphere was 
about 2.4 instead of 3.0 and for a long cylinder oriented with axis normal to the field 
was about 1 .5 instead of 2.0. These measurements were made as close as possible 
to the known analytic locations of maximum field enhancement. At this time, a copper 
coated B-B shot (diameter .173") was being utilized for the sensing probe. 

In order to gain some insight into the reason for the errors in measurement, 
tests were done with varying sizes of probe spheres against a fixed diameter test 
sphere of 3.3 inches. Also varying diameters of both spheres and cylinders were 
sampled with a fixed diameter probe sphere. The results are given in Tables 2.1 - 2.3. 
The measurements on the spheres and cylinders are shown plotted in Figure 2.10. 

The test results indicate that, at least for the large sphere, there is only a 
very small dependence of the measurements on the probe diameter. However, it is 
seen that there is considerable variation of the measurements with the diameters of 
the test spheres and cylinders, especially as the diameters get smaller. 
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Table 2.1 

Form Factors for Fixed Sphere with Varying Probe 
Diameter - Test Sphere Diameter Is 3.3” 


Probe Diameter 

Measured Form Factor 

.092" 

2.68 

.125" 

2.68 

.173" 

2.64 

.30" 

2.63 


Table 2.2 

Form Factors for Varying Sphere Diameters But Fixed 

Probe (.125" dia.) 


Sphere Diameter 

Measured Form Factor 

3.3" 

2.68 

2.0" 

2.53 

1.5" 

2.44 

1.2" 

2.38 


Table 2.3 

Form Factors for Varying Cylinder Diameters but Fixed 

Probe (.125" dia.) 

Long Cylinders Normal to Field 


Cylinder Diameter 

Measured Form Factor 

3.1" 

1.63 

2.5" 

1.59 

1.6" 

1.54 

1.1" 

1.45 

.88" 

1.40 
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Figure 2.10 Correlation Between Measurements and Analytical and Finite 
Difference Solutions for Form Factors 
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Apparently, the form factor measurements move closer to the theoretical values as the 
radius of curvature of the test article at the contact point increases. 

In order to check the observations more carefully, measurements were made 
on a limited size flat plate (8" x 8") and a 2.5 inch cube. The flat plate was oriented 
with its surface normal to the field direction, and the cube with one of its coordinate 
axes along the field direction. The measurement of the flat plate was made in the 
center of the face, and for the cube in the center of one of the faces perpendicular to 
the ambient field. The flat plate measured 0.97 (only a 3% error) and the cube 
measured 1 .63 at the center of the faces normal to the field. A finite difference 
computation of the form factor of a cube in a uniform field gave a value of 1 .7 which 
indicates an error of about 4% for the measurement. 

Further checks of the measurement technique were made utilizing a finite 
length cylinder (15" long) oriented along the field. This configuration provides a range 
of form factors from 0 at the center of the cylinder to in excess of 20 at the ends. The 
cylinder problem was first solved by finite difference methods utilizing a very fine grid 
(~ 1/32”). These results are shown in the first column of Table 2.4. The second 
column of the table shows the same solution utilizing a 1/4" grid. Both these cases 
were run with an idealized uniform ambient field. Note the large degradation of the 
computed values when a larger grid spacing is used. The third column in the table 
shows the results of the finite difference analysis where the test geometry is introduced 
(i.e., the end plates and the guard rings). Obviously the plates and guard rings are not 
reproducing the uniform field exactly. However, the deviations are not large and are 
not the main source of the error. 
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Table 2.4 

Computations and Measurements of the Form Factors Along the Axis of a 
Finite Cylinder Oriented Parallel to the Electric Field 




! 



Finite Difference Solutions 

Plates & Guard Rings 


Distance from 



Measurements 

End (inches) 

1/32" grid 

1/4" grid 

In place - 1/4” Grid 

.125” prob dia. 

.25 

12.79 

8.17 

7.36 

7.1 

.50 

11.11 

6.44 

5.79 

5.54 

1.0 

9.75 

4.9 

4.39 

4.23 

2.0 

8.13 

3.55 

3.16 

3.06 

4.0 

5.17 

2.06 

1.80 

1.80 

6.0 

2.22 

0.89 

0.77 

0.74 

7.5 

0 

0 

0 

0 


The last column in Table 2.4 shows the measured values of form factors 
along the cylinder using the .125 inch diameter probe. Note that these measurements 
agree quite well with computations made by modeling the measurement geometry 
with a 1/4" grid (a smaller grid was not possible for this model because of the small 
time step involved and the resulting increase in computing time). 

The results of the analysis and measurements of the finite cylinder indicate 
clearly that the measurement technique is far from accurate. Errors range from 45% 
near the ends to 300% near the center of the cylinder. A correlation plot of all the 
computed and measured form factors obtained during the study are shown in Figure 
2.10. Note that the correlation curve for the cylinder oriented along the field is 
considerably different than for the spheres and cylinders oriented normal to the field. It 
was originally thought that a universal correlation curve could be found that would 
permit translating measured form factor values to actual values. However, it is seen 
from Figure 2.10 that the uncertainty factor is large, even for the simple geometries 
tested. In other words, the charge transferred to a small metal probe touched to a 
conductive object immersed in an ambient field depends not only on the shape and 
size of the probe and the magnitude of the ambient field but it also depends in a 
complex way on the local distribution of the field at the contact points. 
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2.5 Some Theoretical Support for the Measurement Anomalies 

Hara and Akayaki [7] have solved the problem of the small sphere in contact 
with a flat plate in a uniform field. The charge transfer to the sphere is given by: 

Aq = k£oEr2 (2.1) 

where r = radius of small sphere 
E = ambient uniform field 
k = constant. 

Note that the charge transferred is a function of the radius of the small sphere and the 
field only. This is the condition required for the experimental method being 
investigated to be valid. However, another interesting problem has been analyzed 
that shows the experimental method to be invalid. This is the problem of two spheres 
brought into contact in a uniform field [8]. Here the charge transferred to a small 
sphere from a larger sphere is given by 

Aq = 7EoEr2 (2.2) 

where r = radius of small sphere 
E = ambient electric field 

y = parameter which is a function of r/R, where R is the radius 
of the large sphere. 

The parameter y was computed for a number of values of r/R. It was shown that when 
r/R approached zero, y approached a constant as in equation (2.1 ). Thus the two 
solutions agree in the limiting case of the sphere and the flat plate. 

The analytical results for the case of spheres clearly shows that the 
experimental method is not valid for the simple case of a sphere, since the charge 
transferred will vary with the radius of the sphere being measured whereas the form 
factor of a sphere is independent of its radius. It follows that if the experimental 
method is not valid for spheres, it is not likely to be valid for complex shapes with large 
aspect ratios and multiple curvatures. 
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2.6 


Measurements on the F106B Scale Model 


While it was known that the experimental method being utilized was not 
universally valid, measurements were made on the model aircraft anyway in order to 
obtain data for comparison with the earlier estimates of the form factors generated by 
finite difference analysis. Fourteen different points on the aircraft were measured 
including the three positions where the field meters are located on the real aircraft. 
The locations of the test points are shown in Figures 2.1 1 , 2.12 and 2.13. 

The results of the measurements and some rough estimates of adjusted 
values are given in Table 2.5. The adjusted values can only be very approximate 
since they were obtained by estimating the test point locations as they compared with 
the data given in Figure 2.10. These numbers may not be any more accurate than the 
original measurements. They are included to give the reader some feel for the likely 
errors in the measurements. The three points representing the actual locations of the 
field meters are identified by an asterisk. 


Table 2.5 

Form Factor Measurements and Corrected Estimates 
for Points on the F106B 


Aircraft 
Test Points 

Measurements 
K x K y K z 

Adjusted 

K x 

Measurements 
Ky K z 

* 1 

4.8 

1.4 

0.8 

10.6 

2.6 

1.0 

3 

5.3 

0 

1.9 

11.2 

0 

2.3 

* 4 

3.1 

0 

1.7 

8.0 

0 

2.0 

5 

3.8 

0 

2.4 

9.2 

0 

3.2 

6 

1.3 

0 

2.2 

4.2 

0 

2.8 

7 

0.3 

0.5 

1.2 

1.0 

0.6 

1.4 

9 

0.2 

0.7 

0.7 

0.6 

1.3 

0.8 

10 

2.2 

4.4 

0.8 

6.2 

8.6 

0.9 

11 

0.5 

0.5 

0.4 

1.2 

0.6 

0.4 

12 

6.2 

0.9 

4.8 

12.2 

1.6 

7.6 

* 13 

0.5 

0.2 

1.4 

1.2 

0.3 

1.8 

14 

0.5 

1.0 

1.0 

1.2 

1.8 

1.1 
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Table 2.6 shows a comparison of the computed values obtained by finite 
difference analysis, the actual measured values and the adjusted values. Note that 
the measured values compare quite well with the computed values. It is believed that 
this result is obtained because both methods suffer from a similar problem, i.e., both 
methods produce a value that is averaged over a volume in which the field can vary a 
great deal. As a result, the perceived form factors are different than actual, especially 
where they are much larger than unity. 

The estimated corrections in Table 2.6 are given primarily to indicate that there 
may well be large errors in the values now being used, especially the values for K-| X 
and K 3x . 


Table 2.6 

Comparisons of the Form Factors Obtained by Direct Measurement, 
by Finite Difference Analysis, and by Adjusting the Measured Values 

According to Correlation Data 


Points 

Measurements 

Finite Difference Estimates 

Adjusted Measurements 


K x 

Ky 

K z 

K x 

Ky 

K z 

K x 

Ky 

K z 

Ei 

4.8 

1.4 

0.8 

4.38 

1.46 

0.08 

10.6 

2.6 

1.0 

e 3 

3.11 

0 

1.7 

2.47 

0 

1.52 

8.0 

0 

2.0 

e 4 

0.5 

0.2 

1.4 

0.89 

0.15 

1.03 

1.2 

0.3 

1.8 


2.7 Conclusions 

It has been found that the experimental method for measuring form factors 
utilizing a small scale model immersed in a field and a tiny metal probe has serious 
difficulties. The problems apparently stem from the similarity between the scale of the 
field variations near the skin of the model and the size of the small metal probe. 
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Indications are that if the aircraft model were made very large and the probe kept very 
small, the accuracy of the measurements might be acceptable. However, this action 
would defeat the purpose of the technique which was to facilitate the calibration of field 
meters in a small laboratory environment. 
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CHAPTER 3 


LIGHTNING CHANNEL MODELING 

3.1 Introduction 

Lightning discharges that attach to the FI 06 aircraft often show a large 
number of pulses occurring in the range of several microseconds to tens of 
milliseconds or even hundreds of milliseconds following the initial attachment. The 
number of pulses following the initial discharge varies considerably from flash to flash. 
A crude visual count of the FI 06 analog data for the year 1984 shows an average 
repetition rate of about 250 pulses per second. This rate is undoubtedly the lower limit 
to the true rate since the data used are mostly in the form of low resolution 25 cm/sec 
records. The better resolution 50 cm/sec records, where available, invariably show a 
larger number of pulses in the same time interval. Indeed a higher rate of 1000 pulses 
per second has been reported for a lightning strike on the CV-580 aircraft [9]. These 
pulses are probably due to the discharging of charge clusters in the cloud by the 
leader of a lengthening lightning channel. This would explain the variable pulse 
repetition rate and the time delay between consecutive pulses. An efficient lightning 
channel model is desirable for a better understanding of this phenomenon. 

3.2 Equations Describing the Channel 

The formation and propagation of a lightning channel is a complex 
phenomenon. Many physical processes are involved in the formation of the channel. 
While there are many theoretical models to describe individual stages of channel 
formation, there is presently no model to predict the whole event starting from the 
inception of corona discharge to the highly conducting channel. 

An electrical corona model has been quite successful in modeling the first 
microsecond of a triggered lightning strikes to aircraft [2,3]. The emphasis in the 
application of the model so far has been on the initiation and early time interaction of 
the corona discharge with the aircraft. As a fully three dimensional model with a rather 
elaborate set of fluid equations the existing corona model is not appropriate for 
channel modeling because of its great demand on computer resources. This chapter 
reports on the modification of the corona model to study the characteristics of a leader 
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channel propagating away from the aircraft. This entails extending the analysis to 
cover times greater than a microsecond and distances greater than a few meters from 
the aircraft. An obvious simplification is to assume azimuthal symmetry so that a two 
dimensional model is applicable. Another area of possible simplification is in the 
way the conductivities of the cells are calculated. For channel modeling it is 
reasonable to assume that the conductivity of air away from the channel is low and that 
its variation in time can be neglected. With this assumption only the conductivities of 
the cells forming the channel need to be calculated. With a two dimensional model the 
lightning channel is predetermined to be along the z-axis of an r-z coordinate system. 
The conductivities of these cells are calculated from the charged particle densities and 
mobilities. The use of the fluid conservation equations for mass, momentum and 
energy in the calculation of air conductivity has been described in previous reports 
[2,3]. In the present model the equations for the conservation of momentum are not 
used. Instead the fluid velocities of the charged particles are calculated as the product 
of the respective mobility and electric field, \/ s = ± p s E:. Here the s subscript 

represents species type, and can be e, +, or In the model only motion along the 
z-axis (channel) is considered. A new feature in the present model is the 
inclusion of the vibrational energy reservoir for the neutrals. This and other additional 
features are further explained in the following. 

The rate equations for the particle densities are: 

V • (n e V e ) = Q + Gn e - <x e n 9 - pn e n + + a d n. 

V • (n + V + ) = Q + Gn e - |3n e n + - yn_n + (3.1) 

V • (n.VJ = a e n e - yn_n + - a d n_ 

* Gn e + a d n - + 2 ^ n - + K n + ‘ a e n e 


dn 

“at + 

!n± + 

at + 

an. 
IT + 


at 
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where 


G = is the avalanche rate [2,3,10], 
a Q = the electron attachment rate [2], 

p = the electron-ion recombiation coefficient [11], 

Y = the ion-ion recombination coefficient [2], 
a d = the detachment rate [11] and 

Q = the ambient ionization rate due to cosmic rays [2]. 

The detachment rate is a strong function of the ionic temperature and 

humidity. It is insignificant for ionic temperature lower than 1000° K. For an ionic 
temperature of 2000° K, the detachment lifetime of 0 2 is about 31 nanoseconds. For a 

humidity of 4 gram/m 3 , the detachment lifetime increases rapidly to microseconds for 
secondary types like 0 2 - H 2 0. The other coefficients in the rate equations have been 

discussed in detail in the previous reports. The equation for the neutrals is included to 
handle the situation of a fully ionized gas. However it never needs to be solved. The 
neutral density can be obtained simply from the charged particle densities, the 
conservation of mass and the assumption of singly ionized ions. As mentioned earlier, 
the fluid velocities, V s , of the charged particles are calculated from the mobility, p s , by 
\/ s = ± p s E:, where s stands for species (e,+, -). The mobilities are given by 


m e^en + v e + ) 
2.5 x 10' 4 


(3.2) 


where u en 

is the collision frequency for momentum transfer between 
electrons and neutral molecules [12]; 


is the momentum transfer frequency for electron-ion 
collisions [11], 

Pr 

the relative air density, 

% = 

the electronic charge. 

He. H- = 

mobilities. 
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Through elastic and inelastic collisions, energy is transferred from the 
electrons to the neutral molecules and ions. The energy transferred to the neutrals is 
used in the ionization or dissociation of the molecules, in increasing the translational 
energy or in the excitation of the electronic, rotational or vibrational energy states of 
the molecules. The relaxation of the rotational energy states is fast. Therefore the 
energy released in their relaxation can be combined with the gain in the translational 
energy to increase directly the temperature of the neutrals. In contrast the relaxation 
time of the vibrational states is relatively slow. Therefore a separate conservation 
equation for the vibrational energy reservoir, W v , with a temperature of T v is used in 
the model. The ion and neutral energy equations are combined. This can be done 
because both the ions and neutrals are approximately at the same temperature due to 
strong coupling from collisions. The full set of energy density equations is then: 

+ V Ve e = ' Gn e e ion " (« e + P n +) e e < 3 - 3a ) 


+ e q + a d n T n 





= q+( n + + nj E* • V + + (1 - f v ) C en 


+ C 


e+ 


(3.3b) 


Wy (T v ) ■ w v (T n ) 

T vt 


dW W (T ) - W (T ) 

_ v _ x n v v v l v v r v 

dt " T v u en ' r 


(3.3c) 


In these equations, e jon is the ionization energy of a neutral molecule [2]. 
The temperature for a given species is simply defined by e = nkT, where k is the 
Boltzmann constant. The rate of energy exchange between electrons and neutrals is 

< 34 > 
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where u en is the electron-neutral energy exchange collision frequency [12]. The 

excitation of vibrational energy states of the neutral molecules provides an efficient 
path for the transfer of energy from the electrons to the much heavier neutrals. 
Because of this transfer the temperature of the electon fluid is limited to less than 
25,000° K when the neutral temperature is low compared to that of the electron. The 
rate of energy exchange between electrons and ions through coulomb collisions is 

m 0 

C e + = n e 3 IT ^e+ k ( T e ’ T n> ■ ( 3 -5) 


Thus in the electron energy density equation (Equation 3.3a) the first term on 
the right is energy gained from the ambient electric field. The second term is the 
energy lost in ionizing collisions. The third term is the energy lost in attachment and 
recombination. The fourth term is the energy gained from ambient ionization due to 
cosmic rays. The fifth term is the energy gained from detachment, assuming the freed 
electrons are at the ionic temperature. The last two terms are the energy exchanged 
between electrons and the heavy particles in nonionizing collisions. The energy 
density stored in the vibrational reservoir of the neutrals is 


w v 0" v ) 


e c v 


n n e v 

/kTT 


- 1 


(3.6) 


where e v is the average vibrational energy level. It is taken to be equal to the energy 
of the first excited vibrational state of nitrogen molecule, e v =0.28eV [11]. In the 

equations for the energy balance for the heavy particles and the vibrational energy 
reservoir, f y is the fraction of the energy transfer in electron-neutral collisions that goes 

into the excitation of the vibrational states. For a reduced field, E/n n , of ! ess than 

-19 2 

5x10 VM , f v is larger than 0.75 [11]. Therefore for lightning strike environments 
most of the energy goes into the excitation of the vibrational states, is the 

vibrational relaxation time. It is a function of the neutral temperature and humidity. For 
a neutral temperature of 1000° K, is about 40 microseconds and increases rapidly 

as the temperature drops. It is 200 microseconds for a temperature of 400° K [1 1 ]. 
Thus the transfer of the vibrational energy to the translational energy is rather slow at 
low temperature. There are other energy reservoirs which are not included in the 
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energy equations. Among them is the radiation energy of the photons, which may be 
important in the propagation of a lightning channel via photoionization at the head of 
the leader channel. 

In the solution of the conservation equations, the functional dependence and 
the values of the rate coefficients are taken from the references cited. They are in the 
form of curve fitted expressions, graphs or tables. These data are usually measured 
in conditions which are not those found in a lightning environment. Nevertheless for 
lack of better data they are assumed to be applicable to the study of air breakdown. 
From the densities and mobilities of the particles the conductivity for each cell in the 
finite difference grid can be calculated from 

= q(n e p e + (n + + n_ ) p.) . (3.7) 


3.3 Results 

The model has been utilized to examine the propagation of a triggered 
lightning leader from an airplane. The problem space dimensions are 400 meters in 
the z-direction with a resolution of 2 meters, and 150 meters in the radial direction with 
a resolution of 3 meters. In the two dimensional model the gridding for the nose of the 
aircraft is a tapered rod along the z-axis. For the results presented, a 12-cell cylinder 
with zero resistance is used to model the aircraft. The relative air density used is 0.5. 
The initial temperature is 300° K. 

For an ambient field of 0.45 MV/m, the conductivity of the leader channel at 
selected distances from the aircraft is shown in Figure 3.1 . The variation of the 
conductivity of the cell 2 meters from the aircraft as a function of time is shown in 
curve A of Figure 3.1 . Similarly curves B,C, and D show the conductivities of the cells 
at distances of 42 m, 82 m and 112 m from the aircraft respectively. It is seen that air 
breakdown in the cell at a distance of 42 m away from the aircraft (curve B) occurs 
more than .6 microseconds later than that near the nose. The occurrence of air 
breakdown 3 microseconds later at a distance of 1 12 m from the aircraft is shown by 
curve D. Thus Figure 3.1 clearly shows the formation and propagation of the leader 
channel away from the aircraft. For a given cell the conductivity increases very rapidly 
at the onset of the breakdown and then decreases gradually as the high breakdown 
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field passes on to the next cell in front of the leader tip. The breakdown field is 
propagated forward due to a highly conducting channel behind the tip. As the 
conductivity of the channel decreases due to the attachment of the electrons the field 
at the tip of the channel decreases until it becomes less than the breakdown field 
strength whereupon the propagation of the channel is terminated. Thus in this model 
the conductivity of the channel determines the propagation of the channel. The 
conductivity of the channel and hence the length of the channel is strongly dependent 
on the ambient field. For this particular run, with an ambient field of 0.45 MV/m, the 
length of the leader channel is about 112 meters. The corresponding channel lengths 
are 1 0 m and 35 m for ambient fields of 0.25 and 0.35 MV/m respectively. For an 
ambient field of 0.5 MV/m the leader channel reaches the problem space boundary 
and presumably would be able to reach the charge centers that produced the ambient 
field. This result is consistent with observations showing maximum measured 
thunderstorm fields of about 0.4 MV/m. Fields greater than this value would not exist 
because of corona growth around thunderstorm particles. 

An interesting feature of the leader channel is the pulses that originate from 
the cell adjacent to the aircraft as shown by curve A in Figure 3.1 . These pulses 
propagate away from the aircraft with diminishing amplitude. It is seen that the rate of 
rise of the pulses is very fast while their decay is comparatively slow. This is 
characteristic of air breakdown and electron attachment. The phenomenon can be 
explained by the value of the local field at the particular cell of interest. After the 
breakdown wave passes through, the conductivity of the cell starts to decrease. As the 
resistance rises, the field in the cell also rises until another breakdown of air occurs in 
the cell, whereupon the process is repeated to give the pulse structure in the 
conductivity curve. However, the pulse structure described is most probably a 
numerical phenomenon arising from the way the conductivity is calculated. In the 
present model, the conductivity of air is calculated only for those cells that constitute 
the lightning leader channel. Furthermore, to model the motion of the charges, which 
is mainly along the channel, the conductivity associated with the z-component of the 
field is computed as a function of time, while the conductivity associated with the 
r-component (transverse to the channel) is not calculated as a function of time. It is 
kept at a low constant value equal to the air conductivity outside the channel. In a 
previous three dimensional study where the conductivity of air is calculated in every 
cell and where the motion of charges is calculated consistently, no such pulse 
structure is seen in the conductivity curve adjacent to the nose of the aircraft. 
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Therefore the physical conductivity curve should be approximated by the smoothed 
average of the numerical curve. 

From Figure 3.1, the speed of the leader channel is calculated to be 3.6x1 0 7 
m/s. The leader speed decreases to 3 xl 0 7 m/s and 1 .8x1 0 7 m/s for ambient fields of 
0.35 and 0.25 MV/m respectively. It is seen that the leader velocity calculated from the 
model is strongly dependent on the ambient fields and is much higher than the 
measured value of about 10 6 m/s. In this model the locally enhanced field adjacent to 
the aircraft is responsible for the initiation of the corona discharge and subsequent 
propagation of the leader. With the existing radial resolution of 3 meters the 
enhancement factor is low. Consequently a relatively large ambient field is required 
for the initiation of the leader channel. With a larger enhancement factor from a finer 
grid, the ambient field required would be smaller, thereby lowering the calculated 
leader velocity. Numerically for the leader to propagate a few meters a field higher 
than the observed value is also required. This is not surprising because in this model 
the conductivity of the channel determines the field and hence the avalanche rate at 
the tip of the channel. In this respect the present model can be improved by the 
inclusion of other physical processes, such as photoionization and streamers, which 
are active at the tip of the channel. This is an area where further research work is 
much desired and can be fruitful. An initial effort to incorporate the effects of 
streamers has been done. It is rather difficult to include streamers into the finite 
difference code because of their microscopic dimension compared to the cell size in 
the code. Nevertheless the effects of streamers have been incorporated heuristically 
as shown in a previous report [4]. As formulated in reference [4], there are two free 
parameters to characterize the streamers: their velocity and their growth factor. The 
product of the growth factor and the excess positive ion density is the rate of electron- 
ion pairs produced over and above that produced by field and temperature ionizations. 
In the present implementation streamers are included in the cells at the leader tip only. 
When streamers are included the length of the leader channel typically increases by a 
factor of two or three over that without streamers. For the case of streamers with a 
velocity of 10 6 m/s and a growth factor of 10 s , the leader channel would reach the 
problem space boundary when the ambient field is 0.45 MV/m. 

The case in which charge centers exist along the leader channel has 
also been examined by placing charges in selected cells along the path of the 
channel, i.e. the axis. For the cell size utilized in the present model the charges that 
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can be placed in a given cell without exceeding the breakdown field is of the order of 
microcoulombs. With this charge at a distance at least 50m away from the aircraft, the 
field at the aircraft is very much smaller than the breakdown field. Additional ambient 
field of the order of breakdown field has to be imposed in the problem space. As such 
the development of the leader channel is then controlled mainly by the imposed 
ambient field. The effects of the charge centers are confined to the breakdown of the 
cells where the charges reside. Thus the present model needs to be improved to take 
into account the role of charge centers in leader propagation. 

Thus far the emphasis has been on the propagation of the channel. It is also 
of interest to examine the other properties associated with the channel. For an 
ambient field of 0.45 MV/m the nose current as a function of time is shown in Figure 
3.2. The maximum current is 9.6 kA. The nose current is taken to be the current 
in the cell adjacent to the aircraft. The current in a cell is calculated by the product of 
the field and the conductivity of the cell. The current in the channel decreases rapidly 
as the initial ambient field decreases. The corresponding maximum nose current 
decreases from 9.6 kA to 6.8 kA to 4.9 kA as the field decreases from .45 to .35 to .25 
MV/m. The rate of rise of the lightning current from the model is about 1 0 12 A/s. 
However the breakdown of air in a given cell typically occurs in a couple of time steps 
in the model. Consequently the rate of rise of the lightning current is not accurately 
calculated by the present model, at least near the initiation point. 

The densities of the charged particles in the channel six meters away from 
the aircraft are shown in Figures 3.3 - 3.5. The densities rise rapidly to a maximum as 
air breakdown commences in the cell upon the arrival of the leader tip. After the 
passage of the leader tip the densities decrease. The rate of decrease is much slower 
than the avalanche rate, and decreases with time. The rate of decrease is greatest for 
the electron density, especially in the short period immediately following air 
breakdown. The loss of the electrons is mainly through attachment to neutrals, the 
attachment rate is so high that a microsecond after air breakdown starts most of the 
negative charges are carried by negative ions and not by electrons. Since the mobility 
of the negative ions is much lower than electrons, the attachment of the electrons to 
neutrals to form negative ions results in rapid reduction of the conductivity of the leader 
channel. The temperature of the electrons rises rapidly and tends to oscillate around 
20,000°K, which is close to the observed value. The temperature of the heavy 
particles increases by only a few degrees above the initial value, while the 
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temperature of the vibrational reservoir has almost doubled. Therefore most of the 
electronic energy has gone into the excitation of the vibrational states. The transfer of 
vibrational energy to translation energy is hindered by the low neutral temperature. 
The low gas temperature is also a result of low degree of ionization, since the energy 
gained from the electric field by the relatively small number of charge carriers has to 
be shared among the gas particles. At this gas temperature the detachment rate is so 
low that most of the electrons remain attached. It is seen from Figures 3.3 -3.5 that 
after the rapid air breakdown the electron density is at least three orders of magnitude 
down from the ion densities. 

In conclusion the corona discharge model has been successfully modified to 
study the propagation of the leader channel. The model is based on the consideration 
of some of the pertinent physical processes occurring in the breakdown of air; such as 
field and temperature ionization; recombination of charged particles to form neutrals, 
attachment of electrons to neutrals to form negative ions, detachment of electrons from 
negative ions, energy transfer between the various species and the vibrational energy 
reservoir. The length of the channel, the speed of the leader and the current have 
been shown to depend strongly on the strength of the ambient field. It was found that 
the energy gained by the heavy particles is stored in the vibrational energy reservoir of 
the neutrals. As a result the temperature of the heavy particles hardly increases. 
Electrons that were attached remain attached because of the low ionic temperature. 
With the loss of the free electrons the conductivity of the leader channel decreases 
rapidly and its propagation is terminated. 
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CHAPTER 4 


1984 DATA ANALYSIS AND INTERNAL WIRE MODELING 


4.1 Introduction 

The 1 984 thunderstorm season yielded a large quantity of measured 
lightning responses on the F106B. This chapter presents the analysis of that data and 
concentrates on the derivation of lightning currents from the response waveforms 
using the transfer function technique documented in a previous report [2]. These 
lightning currents are then used to drive the finite difference model of the F106B 
shown in Figure 4.1 . This serves two purposes. The first is to verify that the transfer 
function model chosen is correct; essentially this means that the lightning attach point 
is correctly positioned. The verification is done by comparing calculated and 
measured response waveforms for sensors not used in the transfer function process. If 
these compare favorably the model can be trusted, if not, it is likely that the lightning 
attach point must be altered. 

The second purpose of the finite difference model calculation is to determine 
the sources which couple to the inside of the aircraft via seams and apertures. These 
sources can then be used in a separate calculation to determine internal cable 
currents and voltages. The details of how this is done are presented later in this 
chapter. A number of overlays of measured versus predicted internal wire responses 
are presented as well as predicted responses in cases where no measured response 
was available. 

4.2 Lightning Current Analysis 

There are several steps in the procedure which calculates a lightning current 
from a measured response. The procedure, described in more detail in Reference [2], 
is briefly as follows. A finite difference model of the F106B aircraft and surrounding 
space is constructed as shown in Figure 4.2. A lightning channel is added to the 
model in the form of a thin wire running from an assumed lightning attach point 
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Figure 4.1 Block Model of the F106B Indicating Approximate Sensor Locations 
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Figure 4.2 Typical Model Geometry for Triggered Lightning 
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(usually the nose) to the boundary of the problem space. A source at the attach point 
of the thin wire drives a current on the wire simulating the effect of a triggered strike on 
the aircraft. The resulting current injected on the aircraft is Fourier transformed and 
divided into the Fourier transform of each of the calculated sensor responses on the 
F106B. This gives a number of transfer functions relating injected current at a given 
attach point to sensor response. 

The next step is to Fourier transform a measured sensor response from an 
actual lightning strike on the F106B and to divide it by the appropriate transfer function 
calculated earlier. This results in a frequency domain function which is the injected 
current necessary to cause the measured sensor response in the finite difference 
model. This function is inverse transformed to give the time domain current at the 
attach point of the lightning channel. This derived current is then used in the finite 
difference model to calculate both sensor responses and the fields necessary to derive 
sources for internal coupling purposes. Clearly the measured sensor response must 
by definition match the one calculated by the finite difference code, but the other 
sensor responses are not required to do so. If these comparisons are also good, it is 
an indication that the model used is an accurate simulation of the triggered lightning 
event. 


Sixty-five lightning strikes from the 1984 data set were analyzed in this way. 
In all cases the measured sensor response used to derive the lightning current was 
from the B-dot longitudinal or the D-dot forward sensor. The choice between these two 
was based first of all on which was available, if both were not, and secondly on which 
measured waveform gave better results. 

In all cases a predicted quantity for the l-dot sensor response was calculated 
from the lightning current. This was calculated as shown below: 

Al 

l-dot nose = — (4.1) 

At 

where I is the lightning current at the nose sensor, and At is the time step used in the 
finite difference code (1 nanosecond). Predicted responses for l-dot nose 
weregenerated for each 1984 strike on which the sensors functioned properly. In 
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Figures 4.3(a) - 4.3(k) are shown overlaid measured and predicted l-dot nose 
waveforms. Each of these cases was modeled to 1 microsecond. Figure 4.4 shows 
l-dot nose predictions for possible multiple strike events. In these cases the modeled 
responses were calculated to 2 microseconds. 

In summary it may be said that the comparison of the predicted versus the 
measured l-dot nose responses was in general quite good. Waveforms tend to 
overlay in terms of times of minima and maxima, and amplitudes are usually within 
twenty percent. There were several cases in which disagreements were as much as a 
factor of three between measured and predicted values. In these cases it is likely that 
the attach point chosen for the model was incorrect. 

Some general observations about the calculated and measured sensor 
responses may be made. Magnitudes for the l-dot nose sensor responses usually 
ranged ± 25 x 1 0 9 A/s, but were measured as high as ± 50 x 1 0 9 A/s for a few cases, 
l-dot nose response magnitudes of ± 30 x 10 9 A/s in general correspond to D-dot 
forward and B-dot longitudinal magnitudes of ± 15 A/m and ± 1500 T/s, respectively. 
The relationship between these different sensor responses seems to be linear, so it 
can be expected that every 10 10 A/s increase in l-dot nose will result in increases in 
D-dot forward and B-dot longitudinal of 5 A/m and 500 T/s, respectively. It should be 
noted that the maximum values discussed here are those derived from transient 
measurements. Peak detector measurements in most cases found significantly higher 
values. 


4.3 Internal Coupling and Internal Wire Modeling 

Electromagnetic energy generated by lightning channels attached to an 
aircraft can produce fields inside the aircraft. External fields on the aircraft surface can 
couple to the interior via external wiring and antennas leading to the interior, diffusion 
through the aircraft surface, seams and joints, and apertures. One of the most 
significant consequences of internal coupling is that the external fields caused by the 
attached channel can induce voltages and currents on internal wiring routed across or 
near these points of entry into the aircraft. For the FI 06B, the largest mode of internal 
coupling to the fuselage and wing wires appears to be through seams around the 
wheel well and instrument bay doors. 

(text continued on page 62) 
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RMPS/SEC 


TIME (1 X 10E-6 Seconds) 

I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-024, run 001, strike 001 


TIME (1 X 10E-6 Seconds) 

I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-024, run 003, strike 003 


TIME (1 X 10E-6 Seconds) 

I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-025, run 001, strike 001 


TIME (1 X 10E-6 Seconds] 

1-DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-025, run 002, strike 003 


Figure 4.3(a) Measured versus Predicted l-dot Nose Waveforms 







RMPS/SEC 



I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-025, run 003, strike 004 



I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-025, run 004, strike 005 



I -DOT NOSE 

- MEASURED; -- PREDICTED 
Flight 84-025, run 007, strike 010 



I -DOT NOSE 

- MEASURED: - PREDICTED 
Flight 84-025, run 008, strike OH 


Figure 4.3(b) Measured versus Predicted l-dot Nose Waveforms 
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RMP5/SEC RMPS/SEC 




I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-025, run Oil, strike 014 


I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-025, run 013, strike 017 



I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-025, run 015, strike 019 



I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-027A, run 001, strike 001 


Figure 4.3(c) Measured versus Predicted l-dot Nose Waveforms 
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AMPS/SEC 


TIME (1 X 10E-6 Seconds) 

I -DOT NOSE 

- MEASUREO; - PREOICTEO 
Flight 84-0278, run 001, strike 001 


TIME (1 X 10E-6 Seconds) 

1-00T NOSE 

- MEASUREO; -- PREDICTED 
Flight 84-027A, run 002, strike 002 



i III i n 


TIME (1 X 10E-6 Seconds) 

I -DOT NOSE 

- MEASURED; -- PREOICTEO 
Flight 84-0278, run 002, strike 002 


TIME (1 X 10E-6 Seconds) 

I-OOT NOSE 

- MEASURED; -- PREDICTED 
Flight 84-027, run 004, strike 005 


Figure 4.3(d) Measured versus Predicted l-dot Nose Waveforms 
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0 .2 .4 .6 .8 1.0 

TIME (1 X 10E-6 Seconds) 

I -DOT NOSE 

- MEASURED; -- PREDICTED 
Flight 84-032, run 004, strike 013 
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0 .2 .4 .6 .8 1 

TIME (1 X 10E-6 Seconds) 

1-DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-035, run 001, strike 001 


Figure 4.3(f) Measured versus Predicted l-dot Nose Waveforms 






RMP5/5EC RMPS/SEC 




TIME (1 X 10E-6 Seconds) 


I -DOT NOSE 

- MEASURED; - PREOICTEO 
Flight 84-035, run 003, strike 004 


I-DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-036, run 001, strike 001 




I-DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-036, run 002, strike 002 


I-DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-036, run 006, strike 009 


Figure 4.3(g) Measured versus Predicted l-dot Nose Waveforms 





AMPS/SEC AMPS/SEC 


TIME (1 X 10E-6 Seconds) 

I -DOT NOSE 

- MEASURED; ~ PREDICTED 
Flight 84-037, run 002, strike 018 


TIME (1 X IOE-6 Seconds) 

I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-037, run 003, strike 022 


TIME (I X 10E-6 Seconds) 

I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-037, run 004, strike 024 


TIME [1 X 10E-6 Seconds) 

I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-037, run 005, strike 029 


Figure 4.3(h) Measured versus Predicted l-dot Nose Waveforms 
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0 .2 .4 .6 .8 1. 

TIME (1 X 10E-6 Seconds) 

I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-037, run 006, strike 034 



0 .2 .4 .6 .8 1. 

TIME (1 X 10E-6 Seconds) 

I -DOT NOSE 

- MEASURED; -- PREDICTED 
Flight 84-037, run 008, strike 041 



0 .2 .4 .6 .8 1. 

TIME (1 X 10E-6 Seconds) 

I -DOT NOSE 

- MEASURED; -- PREDICTED 
Flight 84-037, run 010, strike 052 
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TIME (1 X 10E-6 Seconds) 

I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-037, run 014, strike 058 


Figure 4.3(i) Measured versus Predicted l-dot Nose Waveforms 







AMPS/SEC AMPS/SEC 


TIME (1 X 10E-6 Seconds) 

I -DOT NOSE 

- MEASURED; ~ PREDICTED 
Flight 84-037, run 015, strike 062 


TIME (1 X 10E-6 Seconds) 

I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-037, run 017, strike 066 




0 .2 .4 .6 .8 1.0 

TIME (1 X 10E-6 Seconds! 

I-OOT NOSE 

- MEASURED; -- PREDICTED 
Flight 84-037, run 019, strike 071 
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0 .2 .4 .6 .8 1 

TIME II X 10E-6 Seconds) 

I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-037, run 001, strike 002 


Figure 4.3(j) Measured versus Predicted l-dot Nose Waveforms 








RMPS/SEC RMPS/SEC 




I -DOT NOSE 

- MEASURED; -- PREDICTED 
Flighi 84-042, run 001, strike 001 


I -DOT NOSE 

- MEASURED; - PREDICTED 
Flighi 84-042, run 003, strike 003 




I -DOT NOSE 

- MEASURED; ~ PREDICTED 
Flight 84-042, run 005, strike 005 


I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-043, run 002, strike 003 


Figure 4.3(k) Measured versus Predicted l-dot Nose Waveforms 
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RMPS/SEC 



-lOxlO 9 f 


-15*10- 


.5 1.0 1.5 

TIME (1 X 10E-6 Seconds) 


2.0 


.5 1.0 1.5 


TIME (1 X 10E-6 Seconds) 


I -DOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-037, nun 009, strike 044 


I-OOT NOSE 

- MEASURED; - PREDICTED 
Flight 84-037, run 016, strike 063 



TIME 11 X 10E-6 Seconds) 
PREDICTED I -DOT NOSE 

Flight 84-050, run 07A, strike 16A 


Figure 4.4 l-dot Nose Waveforms for Possible Multiple Strike Events 






Because these seams have no deliberate electrical connections across 
them, it is likely that the electromagnetic coupling through these seams is inductive. 
That is, the internal voltage sources are directly proportional to the time derivative 
ofthe current density on the outside skin of the aircraft. This voltage can be determined 
by calculating or assuming a value for the inductive impedance of the seam. Formally 
it is given by, 


V = L$ (4.2) 

where V is the internal voltage generated at the seam, L is the inductive coupling 
coefficient for the seam, and j is the external current density across the seam. 

For each of the sixty-five cases analyzed in Section 4.2, an internal coupling 
calculation was also done. Fuselage and wing wires were modeled and their voltage 
responses compared to the actual measured responses. Wire voltages were also 
predicted for those lightning events for which wire responses were not recorded or the 
sensors did not function properly. The arrangement of these particular wires with 
respect to the F106B is shown in Figure 4.5. The wires are described as insulated, 
22-gauge wires. The wing wire (test point #100) is grounded near the left wing tip and 
routed inside the leading edge of the wing to the instrument bay. The fuselage wire 
(test point #101) is grounded to the radome bulkhead and routed along the right inside 
of the fuselage to the instrument bay. Both wires are terminated in 50 ohms. 

The wire responses are calculated by modeling the internal wire as a 
transmission line illuminated at various points along its length by electric field sources. 
These sources are calculated from the voltage driver of equation (4.2), which in turn is 
derived from external currents calculated by the F106B finite difference code and an 
assumption about the inductive coupling coefficient of the seams. The transmission 
line model of the internal wire allows the wire to have varying inductance and 
capacitance per unit length. In practice, however, because the wires of interest were 
in large wire bundles and the detailed routing was difficult to follow, a single 
capacitance per unit length was used at all points on the wire. Table 4.1 shows the 
parameters used in the two wire models. 
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1 00 COP-1 Instrumentation From weapons bay along 

Package leading edge of left wing 



Figure 4.5 Fuselage and Wing Wire Configuration in the F106B 


Table 4.1 

Parameters Used in the Transmission Line Models of Internal 

F106B Wires 


Parameter 

Fuselage Wire 

Wing Wire 

Length 

9 meters 

13.5 meters 

Cell Size 

.5 meters 

.5 meters 

External Field 
Source Coupling 
Locations 

Front and back seams 
at wheel well door; 
Front seam at 
instrument bay door 

Front and left side 
seams at the instrument 
bay door 1 

Best Fit Averaged 
Inductive Transfer 
Impedance (L) 2 * * 

7.7 nH-m (for all 3 
seams) 

2.0 nH-m (for both 
seams) 

Propagation Velocity 
for Insulated Wire 

2 x 1 0 8 m/sec 

2 x 1 0 8 m/sec 

Wire Impedance 

50 ohms 

50 ohms 

Wire Capacitance 
per Unit Length 

1 x 10' 10 f/m 

1 x 1 0' 10 f/m 

Wire Inductance 
per Unit Length 

2.5 x 10* 7 H/m 

2.5 x 10‘ 7 H/m 

Terminations 

Shorted at forward 
bulkhead; terminated 
in 50 ohms in the 
instrument bay 

Shorted at left wing 
tip; terminated in 50 
ohms in the instrument 
bay 


Notes: 

1 . The wing wire was assumed to cross the front and left side seams of instrument bay door at 40 and 
50 degree angles, respectively. 

2. The inductive transfer impedances (L) of these seams are unknown. A best fit value was 

determined and used for all the wire computer runs. It is likely that in reality the impedance of these 

seams changes from flight to flight as the doors are opened and closed 


64 











In the figures to follow, each internal wire response was calculated for a 
triggered lightning channel attached at the nose of the F106B and an exit channel at 
the tail. Figures 4.6(a) - 4.6(f) show the internal wire voltage across the 50 ohm 
termination for the cases in which both a measured and predicted waveform are 
available. Figures 4.7(a) - 4.7(t) show the predicted wire voltages for the cases in 
which a measured response was unavailable. Figure 4.8 shows the predicted wire 
voltages for a possible multiple strike event. In this case, the model was run to 2 
microseconds. 

Some general comments may be made about the comparisons between 
measured and predicted internal wire voltages. Many of the waveshapes match quite 
closely, while magnitudes tend to vary more, sometimes by as much as a factor of two. 
However, it is important to note that an averaged value for the inductive transfer 
impedance, L, was used for all fuselage wire models, and another averaged value for 
all wing wire models. It is likely that L can change from flight to flight, as the 
doors are opened and closed. A factor of two variation in magnitude then is not 
surprising. This is supported by the observation that the predicted responses vary 
from the measured responses by approximately the same percentage for all strikes 
occurring on a given flight. 

The largest waveform differences generally appear in the response of the 
wing wire. In particular there are higher frequencies in the predicted response than in 
the measured. These frequencies are related to the length of the wire and correspond 
to reflections, occurring probably on the shorted end of the wire. The fact that these do 
not appear in the measured waveforms may be because of the mutual effects of the 
cable bundles in which the test wire is located. The behavior of a single wire within a 
large bundle of variously sized and terminated wires may be significantly different than 
the behavior of an isolated wire, particularly at the higher frequencies. 

Based on the 1984 flight data, the fuselage and wing wire voltages vary in 
the range ± 35 V and ± 5 V, respectively, as upper limits. Most of the fuselage wire 
voltages are in the range ± 15 V, while the majority of wing wire voltages are in the 
range ±1.5 V. These latter numbers correspond to typical B-dot longitudinal values of 
± 1500 T/s and D-dot forward values of ± 15 A/square meter. Typical peak current 
magnitudes are in the range of ±1000 amperes, but a few strikes were seen to have 
peaks as high as ± 2000 amperes or larger. 
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(text continued on page 93) 



WIRE VOLTRGE (V) WIRE VOLTAGE (V) 




FUSELAGE WIRE VOLTAGES 
- MEASURED; - PREDICTED 
Flight 84-014, run 002, strike 002 


FUSELAGE HIRE VOLTAGES 
- MEASURED; - PREDICTED 
Flight 84-015, run 001, strike 001 




FUSELAGE WIRE VOLTAGES 
- MEASURED; - PREDICTED 
Flight 84-015, run 002, strike 002 


FUSELAGE HIRE VOLTAGES 
- MEASURED; - PREDICTED 
Flight 84-015, run 003, strike 001 


Figure 4.6(a) Measured versus Predicted Internal Wire Responses 


66 





WIRE VOLTRGE (V) WIRE VOLTAGE (V) 



TIME (1 X 10E-S Seconds) 

FUSELAGE WIRE VOLTAGES 
- MEASURED; -- PREDICTED 
Flight 84-015, run 004, neorby 002 


(Note: This labeled nearby event is now 
suspected to have been a direct strike) 



FUSELAGE WIRE VOLTAGES 
- MEASURED; -- PREDICTED 
Flight 84-027, run 002, strike 002 



FUSELAGE WIRE VOLTAGES 
- MEASURED; - PREDICTED 
Flight 84-027, run 001, strike 001 



FUSELAGE WIRE VOLTAGES 
- MEASURED; - PREDICTED 
Flight 84-027, run 004, strike 005 


Figure 4.6(b) Measured versus Predicted Internal Wire Responses 
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1 1 r— t 1 1 1 

m 1 


i j L*| 


r j| 

1 1 1 \ ! 



i J J • 

Mi 1 "' \ 

i V 

*! '! V 

l' v 1 / 

III 

1 If 

i i 

V \ i 

1 Vi ■ 

JJ 

II 

!;» 

II 

II 

1 

i 

i i i i i i i 

1 

[l 

i ■ 

J 1 


TIME (I X 10E-6 Seconds) 

FUSELAGE WIRE VOLTAGES 
- MEASURED; - PREDICTED 
Flight 84-035, run 001, strike 001 


TIME (1 X 1QE-6 Seconds) 

FUSELAGE WIRE VOLTAGES 
- MEASURED; - PREDICTED 
Flight 84-037, run 008, strike 041 


TIME (1 X 10E-6 Seconds) 

TIME (1 X 10E-6 Seconds) 

FUSELAGE WIRE VOLTAGES 

FUSELAGE WIRE VOLTAGES 

- MEASURED; - PREDICTED 

- MEASURED; - PREDICTED 

Flight 84-042, run 001, strike 001 

Flight 84-042, run 003, strike 003 

Figure 4.6(c) Measured versus Predicted Internal Wire Responses 






WIRE VOLTAGE (V) WIRE VOLTAGE (V) 




FUSELAGE WIRE VOLTAGES 
- MEASURED; -- PREOICTEO 
Flight 84-027, run 005, strike 006 


FUSELAGE WIRE VOLTAGES 
- MEASURED; - PREDICTED 
Flight 84-029, run 001, strike 001 




FUSELAGE WIRE VOLTAGES 
- MEASURED; -- PREDICTED 
Flight 84-032, run 002, strike 003 


FUSELAGE WIRE VOLTAGES 
- MEASURED; - PREDICTED 
Flight 84-032, run 003, strike 005 


Figure 4.6(d) Measured versus Predicted Internal Wire Responses 
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WIRE VGLTRGE (V) WIRE VOLTRGE (V) 


TIME (1 X 10E-6 Seconds) 

FUSELAGE WIRE VOLTAGES 
- MEASURED; ~ PREDICTED 
FlighL 84-0^32, run 005, sLrike 005 


TIME (1 X 10E-6 Seconds) 

FUSELAGE HIRE VOLTAGES 
- MEASURED; - PREDICTED 
FlighL 83-033, run 002, strike 003 


TIME (1 X 10E-6 Seconds) 

FUSELAGE WIRE VOLTAGES 
- MEASURED; -- PREDICTED 
FlighL 83-033, run 003, strike 003 


TIME (1 X 10E-6 Seconds) 

FUSELAGE WIRE VOLTAGES 
- MEASURED; - PREDICTED 
Flight 83-050, run 002, strike 002 


Figure 4.6(e) Measured versus Predicted Internal Wire Responses 










WIRE VOLTRGE (V) 


TIME (1 X 10E-6 Seconds) 

WING HIRE VOLTAGES 
- MEASURED; -- PREDICTED 
Flight 84-042, run 001, strike 001 


TIME (1 X 10E-6 Seconds) 

WING HIRE VOLTAGES 
- MEASURED; - PREDICTED 
FlighL 84-042, run 005, strike 005 



0 .2 .4 .6 .8 1.' 

TIME (1 X 10E-6 Seconds) 

WING WIRE VOLTAGES 
- MEASURED; - PREDICTED 
Flight 84-043, run 002, strike 003 



0 .2 .4 .6 .8 1. 

TIME (1 X 10E-6 Seconds) 

WING WIRE VOLTAGES 
- MEASURED; - PREDICTED 
Flight 84-043, run 003, strike 004 


Figure 4.6(f) Measured versus Predicted Internal Wire Responses 


71 







WIRE VOLTRGE (V) 




PREDICTED FUSELAGE HIRE VOLTAGE PREDICTED FUSELAGE WIRE VOLTAGE 




PREDICTED WING WIRE VOLTAGE PREDICTED WING HIRE VOLTAGE 


Flight 81-021, run 001, strike 001 


Flight 81-021, run 003, strike 001 


Figure 4.7(a) Predicted Internal Wire Responses 







WIRE VOLTRGE (V) WIRE VOLTRGE (V) 




PREDICTED FUSELAGE WIRE VOLTAGE PREDICTED FUSELAGE WIRE VOLTAGE 




PREDICTED WING WIRE VOLTAGE PREDICTED WING WIRE VOLTAGE 

Flight 84-021 run 001, strike 001 Flight 84-024, run 003, strike 003 

Figure 4.7(b) Predicted Internal Wire Responses 
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PREDICTED FUSELAGE HIRE VOLTAGE 



PREDICTED WING WIRE VOLTAGE 
Flight 84-025, run 001, strike 001 



PREDICTED FUSELAGE WIRE VOLTAGE 



PREDICTED WING WIRE VOLTAGE 
Flight 84-025, run 002, strike 003 


Figure 4.7(c) Predicted Internal Wire Responses 
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WIRE VOLTAGE (V) WIRE VOLTAGE (V) 



PREDICTED FUSELAGE HIRE VOLTAGE 



PREDICTED WING WIRE VOLTAGE 
Flight 84-025, run 003, strike 004 



PREDICTED FUSELAGE WIRE VOLTAGE 



PREDICTED WING WIRE VOLTAGE 
Flight 84-025, run 004, strike 005 


Figure 4.7(d) Predicted Internal Wire Responses 
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PREDICTED FUSELAGE HIRE VOLTAGE 



PREDICTED HING WIRE VOLTAGE 
Flight 84-025, run 007, strike 010 



PREDICTED FUSELAGE WIRE VOLTAGE 



PREDICTED WING WIRE VOLTAGE 
Flight 84-025, run 008, strike Oil 


Figure 4.7(e) Predicted Internal Wire Responses 
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PREDICTED FUSELAGE HIRE VOLTAGE 



PREDICTED WING HIRE VOLTAGE 
Flight 84-025, run Oil, strike OH 



PREDICTED FUSELAGE HIRE VOLTAGE 



PREDICTED WING WIRE VOLTAGE 
Flight 84-025, run 013, strike 017 


Figure 4.7(f) Predicted Internal Wire Responses 
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WIRE VOLTAGE (V) WIRE VOLTAGE (V) 



PREDICTED FUSELAGE WIRE VOLTAGE 



PREDICTED WING WIRE VOLTAGE 
Flight 84-025, run 015, strike 019 



PREDICTED FUSELAGE WIRE VOLTAGE 



PREDICTED WING WIRE VOLTAGE 
Flight 84-0270, run 001, strike 001 


Figure 4.7(g) Predicted Internal Wire Responses 
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WIRE VOLTRGE (V) WIRE VOLTRGE (V) 




PREDICTED FUSELAGE WIRE VOLTAGE PREDICTED FUSELAGE WIRE VOLTAGE 



PREDICTED WING WIRE VOLTAGE 
Flight 84-032, run 001, strike 002 



PREDICTED WING WIRE VOLTAGE 
Flight 84-032, run 002, strike 004 


Figure 4.7(h) Predicted Internal Wire Responses 
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AGE (V) 










WIRE VO 


0 .2 .6 .8 1.0 

TIME (1 X 10E-6 Seconds) 

PREDICTED WING WIRE VOLTRGE 
Flight 84-036, run 001, strike 001 


-.8 


i i i i . — i — . — i 

0 .2 .4 .6 .8 1.0 


TIME (1 X 10E-6 Seconds) 
PREDICTED WING WIRE VOLTAGE 


Flight 84-036, run 002, strike 002 


Figure 4.7(j) Predicted Internal Wire Responses 
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PREDICTED FUSELAGE HIRE VOLTAGE 



PREDICTED HING HIRE VOLTAGE 
Flight 84-037, run 003, strike 022 



PREDICTED FUSELAGE HIRE VOLTAGE 



PREDICTED HING HIRE VOLTAGE 
Flight 84-037, run 004, strike 024 


Figure 4.7(1) Predicted Internal Wire Responses 
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WIRE VOLTRGE (V) HIRE V OLTRGE (V) 




PREDICTED FUSELAGE WIRE VOLTAGE PREDICTED FUSELAGE WIRE VOLTAGE 




PREDICTED WING WIRE VOLTAGE PREDICTED WING WIRE VOLTAGE 

Flight 84-037, run 005, strike 029 Flight 84-037, run 006, strike 034 


Figure 4.7(m) Predicted Internal Wire Responses 
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WIRE VOLTAGE (V) 












Figure 4.7(o) Predicted Internal Wire Responses 






.0 


TIME (1 X 10E-6 Seconds) 
PREDICTED WING WIRE VOLTAGE 



PREDICTED WING HIRE VOLTAGE 


Flight 84-037, run 015, strike 062 


Flight 84-037, run 017, strike 066 


Figure 4.7(p) Predicted Internal Wire Responses 







VOLTAGE (V) WIRE VOLTAGE (V) 









Figure 4.7(r) Predicted Internal Wire Responses 
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PREDICTED FUSELRGE HIRE VOLTRGE 



PREDICTED WING HIRE VOLTRGE 
Flight 84-050, run 008, strike 020 



PREDICTED FUSELRGE HIRE VOLTRGE • 



PREDICTED HING HIRE VOLTRGE 
Flight 84-050, run 010, strike 025 


Figure 4.7(s) Predicted Internal Wire Responses 
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TIME (1 X 10E-6 Seconds) TIME (1 X 10E-6 Seconds) 

PREDICTED WING WIRE VOLTAGE PREDICTED WING WIRE VOLTAGE 


Flight 84-052, run 001, strike 001 Flight 84-052, run 002, strike 002 


Figure 4.7(t) Predicted Internal Wire Responses 
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»-1 «—t »— I (\J 

TIME (1 X 10E-6 Seconds] 
PREDICTED FUSELAGE WIRE VOLTAGE 



«— * «*-♦ •—« •“* C\1 

TIME (1 X 10E-6 Seconds) 
PREDICTED WING WIRE VOLTAGE 


Flight 84-050, run 07A, strike 16A 


Figure 4.8 Predicted Internal Wire Responses for a Possible 
Multiple Strike Event 
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CHAPTER 5 


THE EFFECT OF THUNDERSTORM PARTICLES ON TRIGGERED 

LIGHTNING 


5.1 Introduction 

The initiation of lightning in a thunderstorm, either triggered or natural, is a 
phenomenon which is not fully understood. Maximum observed electric field levels 
within a cloud are generally significantly lower than the levels needed to start an 
electron avalanche which could lead to a lightning event. It is therefore necessary to 
look for sources of field enhancement which increase the electric field in local areas of 
the thundercloud to values which do result in avalanching. The presence of a large 
conducting object, such as the F106B aircraft, is one way in which this may be 
accomplished. Another candidate for enhanced fields is the region around a naturally 
occurring thunderstorm particle such as a water droplet or ice crystal. These often 
irregularly shaped particles can produce significant enhancements in small 
surrounding volumes. When one considers the number of these particles in a 
thundercloud, it is clear that they can be a major source of electrical corona within the 
cloud, and may contribute to the initiation of lightning. 

In this chapter the effect of precipitation particles on lightning initiation will be 
investigated by studying the small scale alterations that they make in the overall 
electric field distribution of the cloud. It will be assumed that each particle is immersed 
in a uniform field over its volume (a very good approximation), and that individual 
particles do not electrically interact with other particles (possibly a less valid 
assumption). The alteration that the particles cause in the net avalanche rate will be 
calculated, and results from a nonlinear triggered lightning model compared for the 
cases in which thunderstorm particles are and are not present. 

5.2 The Effect of Particles on Corona Formation 

By locally redistributing the electric field energy density in a region of space, 
thunderstorm particles can affect the way in which electrons and ions are produced in 
an electron avalanche. In particular, a uniform electric field which is below breakdown 
intensity may be locally enhanced to a level which produces corona around the 
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particles. This occurs because the avalanche rate is strongly nonlinear with electric 
field, so that small electric field changes can result in large changes in charged 
particle production. 

As an example of this consider the case of a spherical particle in a uniform 
static electric field, E 0 , with magnitude ninety percent of air breakdown intensity. This 
is shown schematically in Figure 5.1. The electrical conductivity of the particle is 
important only for very fast field changes (of the order of microseconds), and in general 
the particle acts as a perfect conductor on the time scale of large scale thunderstorm 
fields. Because of this the well-known analytic development for the altered field 
distribution around a conducting sphere in a uniform electric field may be used. Figure 
5.1 shows the regions around the particle in which fields larger and smaller than the 
ambient field exist. In addition the regions in which the field magnitude exceeds air 
breakdown level, Eb, are shown. These latter regions depend on the magnitude of the 
ambient field with respect to breakdown field. That is, if the ambient field drops, the 
breakdown region will shrink, finally disappearing altogether when the ambient field 
goes below one-third of the breakdown field. The other regions, those of enhanced 
and decreased field, are not dependent on ambient field strength. 

The preceding discussion was based on a spherical particle. The situation 
for more irregular particles such as ice crystals or snow flakes is more complicated. 
Analytic solutions for these shapes do not exist, so approximations must be made in 
specifying the local field distribution around the particle. In general, however, a 
polarized particle will behave like an electric dipole to some level of approximation. 
The assumption used in the development presented here is to treat the local field 
distribution the same as that around a spherical particle, but to raise the maximum field 
enhancement to that determined from a previous study [3] (e.g., enhancement of five 
for columnar ice crystals and nine for a typical snowflake). 

The importance of Figure 5.1 is that it shows that a conducting particle can 
cause electron avalanching to occur in a region of space which would have no 
avalanching in the particle's absence. The details of how much avalanching occurs 
because of the particle and the effects of irregularly shaped particles and particle 
density are left for the mathematical development of the next section. 
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Thunderstorm Particle 

E 0 Is Ambient Field, Oriented Horizontally 

Eb Is the Field Required for Air Breakdown 



5.3 


Mathematical Development 


To facilitate the mathematical development to follow, some symbols will first 
be defined. The development will be done in a spherical coordinate system, assuming 
no azimuthal variation. The ambient field is assumed to be oriented along the polar 
axis. The assumption of no azimuthal variation is strictly true only for spherical 
particles or cylindrical particles oriented with their axis along the polar axis. Particles 
not satisfying this criterion will be included in an approximate sense. 

Table 5.1 

Definition of Symbols Used in Section 5.3 


a 

- approximate half size of particle (radius in the 
case of a spherical particle) 

E o 

- ambient electric field intensity 

E b 

- breakdown electric field intensity 

E r 

- radial electric field intensity 

E e 

- theta component of electric field intensity 

E * 

- component of electric field intensity along polar 
axis 

*max 

- maximum electric field enhancement factor 

(particle dependent, defined as the ratio of local 
electric field to ambient electric field) 

n 

- density of thunderstorm particles 

N 

- number of thunderstorm particles in a finite 
difference cell 

P 

- air density relative to sea level density 

AV 

- volume around a particle in which field is above E b 
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The objective of the analysis to follow is to develop a modified electron 
avalanche rate which accounts for the presence of thunderstorm particles. Previous 
work has assumed that corona formation occurs in clear air only. To include particles, 
it is necessary first to calculate the avalanche occurring around a single particle and 
then to consider the effects of a particle density in a given region of space. The electric 
field around a single particle is investigated first. 

The field distribution outside of a spherical conducting particle placed in a 
uniform electric field is given by, 



E e (r,e) = -E 0 (l--j§-)sin9 


Figure 5.2 shows the spherical coordinate system used in the analysis. For use here 
this expression is generalized to approximate the field distribution around irregularly 
shaped particles as follows, 


E, (r.6) = E 0 [l +(f max ~V~p~ Jcose (5.2) 

Eg (r.6) = - E 0 (l- -jj- ) sin 8 

Note that the only difference between equations (5.1) and (5.2) is that the 2 in the first 
of equations (5.1) has been replaced by the term f max -1 in equation (5.2). This 

reduces to equation (5.1 ) in the case of a sphere (for which f max = 3). The justification 

for using equation (5.2) in the present model is not rigorous. Any conducting particle 
placed in a uniform electric field must polarize, and the resulting field distribution must 
in some sense resemble a dipole field. The assumption is made here, mainly for the 
purpose of making further analytic progress on particles in general, that the distribution 
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Figure 5.2 Spherical Coordinate System Used in Analysis of 
Thunderstorm Particles 



is similar to that around a spherical particle, but with differing maximum field 
enhancements. 

In addition, the electric field used in the calculations to follow is the 
component along the polar axis of the coordinate system. Strictly speaking, the total 
electric field is more appropriate. However, the form of the polar component is much 
easier to work with analytically, and introduces only a small error. The error is small 
because the axial oriented field is a good approximation to the total field at most 
locations. The axial field is easily calculated from equation (5.2) and is shown below. 

E z (r.6) = E 0 [l + (f max cos 2 0 - 1 )-pp ] (5.3) 


This expression involves three physical quantities of interest, the ambient 
electric field, the size of the particle, and its shape (included in f max ). These three 

quantities, along with the breakdown field strength at a particular altitude and the 
particle density, completely determine the altered avalanche rate from the clear air 

expression. The next step in the analysis is to specify the volume around the particle 
over which avalanching occurs. In other words, over what volume is E z larger than 

E b ? The condition for this is expressed in equation (5.4). 

E„[l ♦Cmax 00520 ' 1 )-^- ] ^ E b (5.4) 


There are three possible regimes for equation (5.4). If E 0 is larger than E b , 

then equation (5.4) is satisfied in all of space. This is of course the trivial case in which 
air breakdown would occur in the absence of particles. A second regime is one in 
which equation (5.4) is satisfied nowhere in space. This occurs when f^ is not large 
enough to raise the maximum enhanced field above E b . In this case air breakdown 

is absent even in the presence of particles. The middle regime in which equation (5.4) 
is satisfied in a limited region of space is the most interesting. The conditions 

fmax > jE^and E o < E b define this regime. Here air breakdown occurs when particles 
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are present but not when they are absent. From equation (5.4) the volume over which 
this occurs can be specified as below. 


0 < 0 <cos 






max^o 


(5.5) 


a < r <a 



a Ba 


Equation (5.5) represents a two dimensional projection of a three 
dimensional volume. The full volume is obtained by rotating the area specified in 
equation (5.5) around the polar axis of the particle. Note also that there are two of 
these volumes, one on either end of a symmetrically shaped particle. An expression 
for the volume over which the enhanced electric field is above E can now be derived 
by integrating over the limits of equation (5.5) and multiplying by two to account for 
both ends of the particle. Formally this is written, 


2x A Ba 0 

AV = 2 j j j rsin 0 drd0dcp 


(5.6) 


o o a 


The details of the integration are straightforward, and the final expression for 
the avalanche volume is, 


AV = 


4*a 3 E„f, 


o max 


9 E h -E n 

D 0 


\ o max/ \ o max/ 


The most important thing to notice about equation (5.7) is the way in which 
AV depends on various quantities. A large AV implies that the particle has a 
significant effect on electron avalanching, and a small AV implies very little effect. 
Therefore it is clear that larger particles (large a) are most effective at altering the 
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avalanche rate, since AV varies as the cube of the particle dimension. The shape is 
also involved, but is less important, in the sense that AV is nearly linearly dependent 
on f max . Of course it must be remembered that f max is quite important in determining 

whether there will be any avalanching at all around the particle. 

Equation (5.7) gives the volume around a single particle in which 
avalanching will occur. To implement this concept in a numerical model, it is 
necessary to include the density of particles which is present in the region of interest. 

In a finite difference model the numerical region of interest is always a finite difference 
cell. Hence it is necessary to calculate the total volume which has fields above 
breakdown level in a given cell. The assumption is made that the electric field 
calculated in a finite difference cell is the ambient field in that cell, and that any particle 
density that may be present acts to enhance that field in small volumes around the 
particles. The total volume with field above breakdown strength in a single finite 
difference cell is then, 


AV t = NAV (5.8) 

where N is the number of particles in that cell and is numerically equal to the particle 
number density multiplied by the cell volume. In equation (5.8) the assumption has 
been made that the particles are noninteracting in the sense that none of the individual 
enhanced field volumes overlap. This appears to be a reasonable approximation for 
typical particle densities in thunderstorms. The assumption may be restated by saying 
that the total volume above breakdown must be less than the total cell volume. 

The next step in the analysis is to calculate the number of electron-ion pairs 
produced in AV t per unit time. This is essentially the desired avalanche rate from 

particles. The rigorous method to find this quantity is to integrate the (strongly electric 
field dependent) avalanche rate over AV. This is computationally inefficient in a finite 
difference code, however, so the method used here is somewhat different. The 
alternate method is to define an average electric field in the breakdown volume and to 
simply evaluate the altered avalanche rate for that Electric field. One possible choice 
of this average field is to set it equal to the arithmetic average of the smallest and 
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largest fields in the volume. These fields are E b and f max E 0 , respectively. However, 
because the largest field exists only in a very small subvolume of AV it would be better 
to favor the lower fields in any average. Therefore the geometric mean is used, as 
shown below. 


Average field = ^E b f max E 0 


(5.9) 


The geometric mean is always less than or equal to the arithmetic mean, so this 
should be a more appropriate value to use in the model. 

Finally the avalanche rate from particles may be formally written, 


AV 


^particles ^clear air( avera 9 e v 


3Dcell 


(5.10) 


This form is suitable for inclusion in the finite difference model documented in previous 
work [3]. It is simply an adjustment to the usual (clear air) avalanche rate by the ratio of 
volumes, and with the usual avalanche rate evaluated at an appropriate average field. 

The last matter to be considered is the issue of irregularly shaped particles 
which may not be aligned with their axis along the ambient electric field. Although 
electrostatic forces would tend to align these particles, it is likely that hydrodynamic 
forces from winds and small scale turbulence would overwhelm the electrostatic force. 
Hence the particles are probably randomly oriented with respect to the field direction. 
This is taken into account in the present model by assuming that on the average half of 
the particles are aligned along the field and half are aligned across the field. The 
particles across the field make no contribution to the particle dependent avalanche 
rate. In the model this is accomplished simply by dividing the true particle density by 

two to arrive at an effective density. It should be understood that this is done only for 
the case of nonspherical particles (f max ^3), because for spherical particles alignment 

is unimportant. 
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5.4 


Results 


One of the immediate results which can be calculated from the thunderstorm 
particle model is the effective ambient field which causes air breakdown to occur. This 
should be lower in the presence of thunderstorm particles than in their absence. An 
approximation to the breakdown field is achieved by finding the field level at which the 
avalanche rate becomes larger than the electron attachment rate. Although this is not 
a strict criterion for air breakdown it does serve to illustrate the difference in corona 
growth caused by particles. Three situations were analyzed: 

(1) No particles 

(2) Small particles with high density 

a = 1 mm 

n = 1 x 1 0 5 particles/cubic meter 

f max= 5 

(3) Large particles with low density 

a = 2 cm 

n = 100 particles/cubic meter 

w- 9 - 

The cases were done at ten different relative air densities ranging from .1 to 1 .0. All 
cases were investigated for an assumed water content of 0 percent. This may not be 
appropriate for the center of a thunderstorm, but the results are not greatly changed for 
any realistic water content. Table 5.2 shows the nominal breakdown field for the three 
cases above as a function of relative air density. The labeling of the cases 
corresponds to the numbering given above. 
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Table 5.2 

Nominal Air Breakdown Field as a Function of Air Density 


p 

E (MV/m) 

Case 1 

Case 2 

Case 3 

.1 

.22 

.21 

.16 

.2 

.43 

.42 

.33 

.3 

.64 

.63 

.49 

.4 

.86 

.84 

.65 

.5 

1.07 

1.05 

.81 

.6 

1.29 

1.26 

.97 

.7 

1.50 

1.48 

1.13 

.8 

1.72 

1.69 

1.30 

.9 

1.93 

1.90 

1.46 

1.0 

2.15 

2.11 

1.62 


The interesting thing to notice from Table 5.2 is that the small high density 
particles have almost no effect on the nominal breakdown field. That case is typical of 
small high density ice crystals that one might find in the upper portions of a 
thundercloud. Case 3 represents quite large particles such as snowflakes. Results 
from the model show that it is necessary to have particles this large in order for the 
particles to have any significant effect on the air breakdown field. Small particles with 
low field enhancement factors have no appreciable effect on the overall corona 
initiation process. Although they enhance the electric field in a local region the volume 
over which this occurs is so small that very few avalanche electrons are produced in a 
macroscopic sense. The reason for this is the cubic dependence on the particle size 
in the expression for the avalanche volume (equation '5.7)). In summary it appears 
that thunderstorm particles are of great significance in altering the breakdown field 
only in regions where large particles exist. 

Although particles do not in most cases greatly change the minimum air 
breakdown electric field, their presence can significantly change the responses seen 
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on an aircraft. For example, if the F106B is placed in a given uniform electric field of 
sufficient magnitude to trigger a lightning strike, the responses seen when particles are 
present are larger than they are in clear air. This is true only for a fixed uniform field 
which is the same in both cases. It can be explained by noting that the breakdown 
process, particularly the electron avalanche rate, is strongly nonlinear with field 
intensity, so that even a small increase of the field strength above minimum 
breakdown can result in a large change in breakdown character. Of course, in the 
case under consideration here, the particles produce a small decrease in the 
breakdown field, equivalent to a small increase in the ambient field. 

The thunderstorm particle model has been implemented in the nonlinear 
finite difference triggered lightning model, and a few different thunderstorm situations 
investigated. One typical of what the F106B may encounter is shown in Figures 5.3 
and 5.4. The case shown there is of a field oriented from top to bottom of the F106B. 
The ambient field magnitude is .3 MV/m, and there is a charge of -.895 millicoulombs 
on the aircraft. The relative air density is .5. Figure 5.3 shows the aircraft responses in 
the case where no particles are present. Figure 5.4 shows the responses for a particle 

environment. The particle density is 1 x 10 5 particles per cubic meter, the particle 
radius is 1 millimeter, and f max is set at 5. This corresponds to columnar ice crystals at 

a high altitude. 

It has been assumed for this comparison that the ambient electric field is the 
same in both the particle and clear air cases. This results in larger aircraft responses 
in the particle case because of the larger avalanche rate. In practice it is likely that the 
presence of particles in the air would slightly lower the triggering field and produce 
aircraft responses similar in magnitude to the pure air case. This would occur because 
the initial lightning current which drives all of the aircraft sensor responses is strongly 
dependent on the temporal behavior of the air conductivity. A slight lowering of the 
ambient electric field slows the growth of the lightning current which in turn drops the 
peak sensor responses seen on the F106B. 

The main things to notice about the comparison are the following. First, the 
magnitudes of the responses for the case with particles are approximately a factor of 
two larger for the D-dot responses and five or more larger for the B-dot response. This 
occurs because of the slight lowering of the effective breakdown field in the 
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Figure 5.3 F106B Responses for a Vertically Oriented Ambient Reid with 

Thunderstorm Particles Absent 
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Figure 5.4 FI 06B Responses for a Vertically Oriented Ambient Field with 
Thunderstorm Particles Present 
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particle case. This is also evidenced by the slightly earlier breakdown time. The 
second thing to notice, which may not be obvious from the plots themselves, is that the 
character of the responses is not greatly altered. If one ignores the larger first peak in 
the particle case, the rest of the waveforms are quite similar in terms of shape. The 
difference in scale of the plots of the two cases tends to obscure this character. The 
similarity is not completely unexpected, as the later portions of the waveforms must be 
characteristic of the resonances of the FI 06B, and are the same whether or not 
particles are present in the surrounding air. Although not shown in the figures, it may 
also be noted that the peak lightning currents in the two cases differed by 
approximately a factor of two, as would be expected from the D-dot responses. 

In summary, the following has been learned from the thunderstorm particle 

model: 


1 . The electric field at which air breakdown occurs is not greatly altered by 
the presence of thunderstorm particles, except for the largest particles. 

2. Particle size is much more important in altering the clear air avalanche 
rate than is particle density. This is because of the cubic variation of the 
avalanche volume with particle size. 

3. Even though the effective breakdown field is not greatly changed by the 
presence of particles, aircraft responses may be significantly affected if 
breakdown occurs at the same field levels. This occurs because of the 
strongly nonlinear variation of breakdown characteristics with electric 
field. The net effect of particles is more likely to be a slight lowering of 
the trigger field with similar aircraft responses. 
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CHAPTER 6 


SUMMARY AND CONCLUSIONS 

The work documented in this report has concentrated on the triggered 
lightning environment of the F106B aircraft and the effect this environment has on 
aircraft internal and external electromagnetic responses. Models have been 
developed for lightning channels and internal wires on the aircraft. The effect of 
particles in the triggering event has been investigated, and form factors of the F106B 
measured. Several significant achievements may be highlighted. 

6.1 Form Factor Measurement 

Although the method used to measure form factors on the F106B was not as 
successful as might be desired for large enhancements, the work performed was 
noteworthy nevertheless. The method has been widely utilized by other investigators 
and cited in the literature. However, no evaluation of the accuracy of the method by 
comparison with canonical shapes has been reported. The comparisons done in the 
present work define the limits within which reasonable results may be expected. Also, 
it has been shown that the method is accurate to 1 0 to 20 percent for moderate field 
enhancements in regions where the electric field is not rapidly varying in space. This 
should be true of the locations of the field mills on the F106B. Hence the form factors 
measured for these field mill locations may be usable in analyzing the measured field 
mill records on the aircraft. However, Mazur [13] has shown that even small errors in 
form factors can translate to large errors in ambient fields so any error in the form 
factors must be treated with caution. 

6.2 Lightning Channel Modeling 

A lightning channel model utilizing a more complete air chemistry has been 
developed. The model, because it is two dimensional, is capable of following the 
development of the channel to late time (tens of microseconds). The model has been 
used to simulate a triggered strike to the nose of the F106B by placing a conical 
conducting structure representing the nose into the problem space. 
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6.3 


Internal Wire Modeling 


The modeling of internal wires on the F106B has been quite successful, 
especially in view of the complicated internal structure of the aircraft. It has been 
shown that simple transmission line models employing discrete sources are adequate 
to simulate the majority of internal wire responses. The success of the modeling also 
provides additional validation of the linear transfer function models used to derive 
sources for the internal wires. These models, which use a measured response 
waveform to derive a lightning current from which internal wire sources are derived, 
have also been successfully used in modeling multiple simultaneous external 
electromagnetic responses. Examples of this are the successful overlays of l-dot 
measured and calculated responses presented in this report. 

6.4 Thunderstorm Particle Modeling 

The modeling done in the investigation of the effect of thunderstorm particles 
on triggered lightning events has shown that the presence of particles does not greatly 
affect the threshold breakdown field intensity. The largest particles were shown to 
have the most effect on the breakdown field. However, because of the strongly 
nonlinear character of the breakdown on electric field strength, even a small change in 
threshold breakdown field can result in a large change in aircraft response. This has 
been illustrated by calculating the responses on the aircraft for both the clear air model 
and the model containing thunderstorm particles. The breakdown was found to be 
significantly larger in terms of aircraft response peak magnitudes for the case having 
particles. This finding could be of value in explaining some of the variability seen in 
measured FI 06B electromagnetic responses. The presence of various types of 
particles having varying shapes and enhancements in the air surrounding the F106B 
introduces an element of randomness into the triggered lightning process that is not 
present for the case of the aircraft flying in clear air. It should also be noted, however, 
that the most likely role for particles to play in triggered lightning is to slightly lower the 
initial breakdown field. The actual electron avalanche should then proceed similarly 
(with respect to the growth rate of the air conductivity), leading to similar aircraft 
responses for the pure air and particle cases in terms of peak amplitudes. 
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